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ABSTRACT 

In order to understand the dynamical and chemical evolution of our Galaxy it is 
of fundamental importance to study the local neighborhood. White dwarf stars are 
ideal candidates to probe the history of the solar neighborhood, since these "fossil" 
stars have very long evolutionary time-scales and, at the same time, their evolution is 
relatively well understood. In fact, the white dwarf luminosity function has been used 
for this purpose by several authors. However, a long standing problem arises from the 
relatively poor statistics of the samples, especially at low luminosities. In this paper 
we assess the statistical reliability of the white dwarf luminosity function by using a 
Monte Carlo approach. 
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1 INTRODUCTION 

The white dwarf luminosity function has become an impor- 
tant tool to determine some properties of the local neigh- 
borhood, such as its age (Winget et al. 1987; Garcia-Berro 
et al. 1988; Hernanz et al. 1994), or the past history of the 
star formation rate (Noh & Scalo 1990; Diaz-Pinto et al. 
1994; Isern et al. 1995a,b). This has been possible because 
now we have improved observational luminosity functions 
(Liebert, Dahn & Monet 1988; Oswalt et al. 1996; Leggett, 
Ruiz & Bergeron 1998) and because we have a better under- 
standing of the physics of white dwarfs and, consequently, 
reliable cooling sequences — at least up to moderately low 
luminosities. 

The most important features of the luminosity function 
of white dwarfs are a smooth increase up to luminosities of 
log(Z//Z/Q) ~ —4.0, and the presence of a pronounced cut- 
off at log(Z//Z/Q) ~ —4.4, although its exact position is still 
today somehow uncertain since it hinges on the statistical 
significance of a small subset of objects, on how the avail- 
able data is binned and on the fine details of the sampling 
procedure. Most of the information on the early times of 
the past history of the local neighborhood is concentrated 
on this uncertain low luminosity portion of the white dwarf 
luminosity function. 

A major drawback of the luminosity function of white 



dwarfs is that it measures the volumetric density of white 
dwarfs and, therefore, in order to compare with the obser- 
vations one must use the volumetric star formation rate, 
that is the star formation rate per cubic parsec, whereas for 
many studies of galactic evolution the star formation rate 
per square parsec is required and, consequently, fitted to 
the observations. 

Another important issue is the fact that the sample 
from which the low luminosity portion (Mv > 13 mag ) of 
the white dwarf luminosity function is derived has been se- 
lected on a kinematical basis (white dwarfs with relatively 
high proper motions). Therefore, some kinematical biases or 
distortions are expected. Although there are some studies of 
the kinematical properties of white dwarf stars — see, for 
instance, Sion et al. (1988) and references therein — a com- 
plete and comprehensive kinematical study of the sample 
used to obtain the white dwarf luminosity function remains 
to be done. It is important to realize that a conventional 
approach to compute theoretical luminosity functions (Her- 
nanz et al. 1994; Wood 1992) does not take into account 
the kinematical properties of the observed sample. A Monte 
Carlo simulation of a model population of white dwarfs is 
expected to allow the biases and effects of sample selection 
to be taken into account, so their luminosity function could 
be corrected — or, at least, correctly interpreted — provided 
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that a detailed simulation from the stage of source selection 
is performed accurately. Of course, a realistic model of the 
evolution of our Galaxy is required for that purpose. 

Finally, the available white dwarf luminosity functions 
(Liebert et al. 1988; Oswalt et al. 1996, Leggett et al. 1988) 
have been obtained using the 1/Vmax method (Schmidt 
1968), which assumes a uniform distribution of the objects, 
yet nothing in our local neighborhood is, strictly speaking, 
homogeneous. In fact, stars in the solar neighborhood are 
concentrated in the plane of the galactic disk. Moreover, it 
is expected that old objects should have larger scale heights 
than young ones. This dependence on the scale height proba- 
bly has effects on the observed white dwarf luminosity func- 
tion — especially on its low luminosity portion where old 
objects concentrate — and, once again, a realistic model of 
galactic evolution is required for evaluating the effects of 
the departures from homogeneity of the observed samples. 
To our knowledge, this effect was only taken into account 
for the bright portion of the white dwarf luminosity func- 
tion (Fleming, Liebert & Green 1986) and not for the low 
luminosity portion where the effects are expected to be more 
dramatic. 

Perhaps the most sucessful application of the white 
dwarf luminosity function has been its invaluable contri- 
bution as an independent galactic chronometer to a better 
understanding of our Galaxy. Despite this fact there have 
been very few attempts — being those of Garcfa-Berro & 
Torres (1997), Wood (1997) and Wood & Oswalt (1998) the 
only serious ones — to systematically investigate the sta- 
tistical uncertainties associated with the derived age of the 
disk. Nevertheless the approach used by Wood & Oswalt 
(1998) makes use of the observed kinematic properties of the 
white dwarf population instead of using a standard model 
of the evolution of our Galaxy, which presumably should 
include the effects of a scale height law. Besides, these au- 
thors use the theoretical white dwarf luminosity function ob- 
tained from standard methods to assign probabilities and, 
ultimately, to assign luminosities to the white dwarfs of the 
sample. Finally, in their calculations Wood & Oswalt (1998) 
computed the cooling times of all the white dwarfs of their 
sample by interpolating in a model cooling sequence of a 
0.6 Mq white dwarf, thus neglecting the effects of the full 
mass spectrum of white dwarfs. 

In this paper we explore the statistical reliability and 
completeness of the white dwarf luminosity function taking 
into account all of the above mentioned effects that were dis- 
regarded in previous studies. Special emphasis will be placed 
on the statistical significance of the reported cut-off in the 
white dwarf luminosity function. For that purpose we will 
use a Monte Carlo method, coupled with bayesian inference 
techniques, within the frame of a consistent model of galac- 
tic evolution, and using improved cooling sequences. 

To be precise, we want specific answers for the follow- 
ing questions: are the kinematics of the derived white dwarf 
population consistent with the observational data? Which 
are the effects of a scale height in the observed samples? Is 
the sample used to derive the white dwarf luminosity func- 
tion representative of the whole white dwarf population? 
Or, at least, is this sample compatible with the white dwarf 
population within the limits imposed by the selection pro- 
cedure? Which are the statistical errors for each luminosity 



bin? Which is the typical sampling error in the derived age 
of the disk? 

The paper is organized as follows: in §2 we describe how 
the simulated population of white dwarfs is built; in §3 we 
describe the kinematical properties of the samples obtained 
in this way and we compare them with those of a real, al- 
though very preliminary and possibly uncomplete, sample; 
in §4 we study the spatial distribution of the samples, we as- 
sess the statistical reliability and completeness of the white 
dwarf luminosity function and we derive an estimate of the 
error budget in the determination of the age of the disk; 
finally, in §5 our results are summarized, followed by con- 
clusions and suggestions for future improvements. 



2 BUILDING THE SAMPLE 

The basic ingredient of any Monte Carlo code is a gener- 
ator of random variables distributed according to a given 
probability density. The simulations described in this paper 
have been done using a random number generator algorithm 
(James 1990) which provides a uniform probability density 
within the interval (0, 1) and ensures a repetition period of 
>, 10 18 , which is virtually infinite for practical simulations. 
When gaussian probability functions are needed we have 
used the Box-Muller algorithm as described in Press et al. 
(1986). 

We randomly choose two numbers for the galactocen- 
tric polar coordinates (r, 6) of each star in the sample within 
approximately 200 pc from the sun, assuming a constant 
surface density. The density changes due to the radial scale 
length of our Galaxy are negligible over the distances we are 
going to consider here and can be completely ignored. Next 
we draw two more pseudo-random numbers: the first for the 
mass (M) on the main sequence of each star — according to 
the initial mass function of Scalo (1998) — and the second 
for the time at which each star was born (t b ) — according 
to a given star formation rate. We have chosen an exponen- 
tially drecreasing star formation rate per unit time and unit 
surface: t/j oc e _t//Ts . This choice of the shape of the star for- 
mation rate is fully consistent with our current understand- 
ing of the chemical evolution of our Galaxy — see, for in- 
stance, Bravo et al. (1993). Once we know the time at which 
each star was born we assign the z coordinate by drawing 
another random number according to an exponential disk 
profile. The scale height of newly formed stars adopted here 
decreases exponentially with time: H p (t) = z\ e~*/ Th + zt. 
This choice for the time dependence of the scale height is 
essentially arbitrary, although it can be considered natu- 
ral. We will however show that using these prescriptions for 
both the surface star formation rate and the scale height law 
to compute the theoretical white dwarf luminosity function 
leads to an excellent fit to the observations (see Isern et al. 
1995a, b and §4) which does not result in a conflict with the 
observed kinematics of the white dwarf population (see §3). 
The values of the free parameters for both the surface star 
formation rate and the scale height have been taken from 
Isern et al. (1995a,b), namely: r s = 24 Gyr, r h = 0.7 Gyr, 
Zl / Z( = 485. 

In order to determine the heliocentric velocities in the 
B3 system, (17, V, W), of each star in the sample three more 
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quantities are drawn according to normal laws: 

n{U) oc e^ u ~ u 'o^'^ 
n(V) oc e - (v - v o' )2/CT v 
n(W) oc e - ( "o )2/CT w 



(1) 



where (Uq, V ', Wq) take into account the differential rotation 
of the disk (Ogorodnikov 1965), and derive from the peculiar 
velocity (Uq, Vq, Wq) of the sun for which we have adopted 
the value (10, 5, 7) km s _1 (Dehnen & Binney 1997). 

The three velocity dispersions (uu, ov, ffw), and the lag 
velocity, Vb, of a given sample of stars are not independent of 
the scale height. From main sequence star counts, Mihalas 
& Binney (1981) obtain the following relations, when the 
velocities are expressed in km s _1 and the scale height is 
expressed in kpc: 
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2 / 2 
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-au/120 



= 



(2) 



0.32 + 1.67 10"Vu 



= 0.50 



6.52 10~Vw 



which is what we adopt here (see as well §3.3). Note, how- 
ever, that our most important input is the scale height law, 
from which most of the kinematical quantities are derived. 

Since white dwarfs are long lived objects the effects of 
the galactic potential on their motion, and therefore on their 
positions and proper motions, can be potentially large, es- 
pecially for very old objects which populate the tail of the 
white dwarf luminosity function. Therefore, the z coordinate 
is integrated using the galactic potential proposed by Flynn 
et al. (1996). This galactic potential includes the contribu- 
tions of the disk, the bulge and the halo, and reproduces 
very well the local disk surface density of matter and the 
rotation curve of our Galaxy. We do not consider the effects 
of the galactic potential in the r and 8 coordinates. This is 
the same as assuming that the number of white dwarfs that 
enter into the sector of the disk that we are considering (the 
local column) is, on average, equal to the number of white 
dwarfs that are leaving it. Of course, with this approach 
we are neglecting the possibility of a global radial flow, and 
thus, the possible effects of diffusion across the disk. How- 
ever, the observed disk kinematics suggest that radial mixing 
is efficient up to distances much larger than the maximum 
distance we have used in our simulations (Carney, Latham 
& Laird 1990). 

From this set of data we can now compute parallaxes 
and proper motions for all the stars (~ 200 000) in the sam- 
ple. Given the age of the disk (tdisk) we can also compute 
how many of these stars have had time to evolve to white 
dwarfs and, given a set of cooling sequences (Salaris et al. 
1997, Garcfa-Berro et al. 1997), what are their luminosities. 
This set of cooling sequences includes the effects of phase 
separation of carbon and oxygen upon crystallization and 
has been computed taking into account detailed chemical 
profiles of the carbon-oxygen binary mixture present in most 
white dwarf interiors. These chemical profiles have been ob- 
tained using the most up to date treatment of the effects 
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Figure 1. Some relevant distributions obtained from a single 
Monte Carlo simulation, the solid lines correspond to the origi- 
nal sample whereas the dotted lines correspond to the restricted 
sample. See text for details. 



of an enhanced reaction rate for the 12 C(a,7) le O reaction. 
Of course, a relationship between the mass on the main se- 
quence and the mass of the resulting white dwarf is needed. 
Main sequence lifetimes must be provided as well. For these 
two relationships we have used those of Iben & Laughlin 
(1989). The size of this new sample of white dwarfs typically 
is of ~ 60 000 stars (hereinafter "original" sample). Finally, 
for all white dwarfs belonging to this sample bolometric cor- 
rections are calculated by interpolating in the atmospheric 
tables of Bergeron et al. (1995) and their V magnitude is 
obtained, assuming that all are non-DA white dwarfs. 

Since the final goal is to compute the white dwarf lu- 
minosity function using the 1/Knax method (Schmidt 1968) 
a set of restrictions is needed for selecting a subset of white 
dwarfs which, in principle, should be representative of the 
whole white dwarf population. We have chosen the follow- 
ing criteria for selecting the final sample: mv < 18.5 mag and 
H > 0.16" yr -1 (Oswalt et al. 1996). We do not consider 
white dwarfs with very small parallaxes (tv < 0.005"), since 
these are unlikely to belong to a realistic observational sam- 
ple. All white dwarfs brighter than My < 13 mag are included 
in the sample, regardless of their proper motions, since the 
the luminosity function of hot white dwarfs has been ob- 
tained from a catalog of spectroscopically identified white 
dwarfs (Green 1980; Fleming et al. 1986) which is assumed 
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to be complete. Additionally all white dwarfs with tangential 
velocities larger than 250 km s" 1 were discarded (Liebert, 
Dahn k, Monet 1989) since these would be probably classi- 
fied as halo members. These restrictions determine the size 
of the the final sample which typically is ~ 200 stars (here- 
inafter "restricted" sample). Finally we normalize the total 
density of white dwarfs obtained in this way to its observed 
value in the solar neighborhood (Oswalt et al. 1996). 

In Figure 1 we show a summary of the most relevant 
results for a disk age of 13 Gyr. In the top panel, the 
mass distribution of those stars that have been able to be- 
come white dwarfs (solid line, left scale) and of those white 
dwarfs that are selected for computing the luminosity func- 
tion (dotted line, right scale) are shown. Both distributions 
are well behaved, follow closely each other, and peak at 
around 0.55 Mq in very good agreement with the obser- 
vations (Bergeron, Saffer & Liebert 1992). In this sense, the 
restricted sample could be considered as representative of 
the whole white dwarf population. 

In the middle panel of figure 1 we show the raw distri- 
bution of luminosities for the stars in the original (solid line, 
left scale) and the restricted (dotted line, right scale) sam- 
ples. The differences between both distributions are quite 
apparent: first, the restricted sample has a broad peak cen- 
tered at log(L/ Lq) ~ —3.5, whereas the original sample 
is narrowly peaked at a smaller luminosity (0.6 dex). Ob- 
viously, since the restricted sample is selected on a kine- 
matical basis — see the lower panel of figure 1, where the 
distribution of proper motions for both samples is shown — 
some very faint and low proper motion white dwarfs are dis- 
carded. Thus, the restricted sample is biased towards larger 
luminosities. Therefore, the cut-off of the observational lumi- 
nosity function should be biased as well towards larger lumi- 
nosities. However, it is important to realize that only ~ 0.6% 
of the total number of white dwarfs with log(L/LQ) > —4.0 
are selected for the restricted sample, and therefore, used in 
computing the white dwarf luminosity function. This num- 
ber decreases to ~ 0.04% if we consider the low luminosity 
portion of the white dwarf luminosity function — that is, 
white dwarfs with log(L/Z/Q) < —4.0 — where most of the 
information regarding the initial phases of our Galaxy is 
recorded. The distribution of proper motions (lower right 
panel of figure 1) shows that most white dwarfs for both 
the original and the restricted sample have proper motions 
smaller than 0.4" yr _1 . However, the restricted sample has 
a pronounced peak at /i ~ 0.3" yr _1 , and shows a deficit 
of very low proper motion white dwarfs, as should be the 
case for a kinematically selected sample, whereas the original 
sample smoothly decreases for increasing proper motions. 



3 THE KINEMATIC PROPERTIES OF THE 
WHITE DWARF POPULATION 

Since the pioneering work of Sion & Liebert (1977), very 
few analysis of the kinematics of the white dwarf popula- 
tion have been done, with that of Sion et al. (1988) being the 
most relevant one, despite the fact that the low luminosity 
portion of the white dwarf luminosity function is actually de- 
rived from a kinematically selected sample. Sion et al. (1988) 
used a specific subset of the proper-motion sample of spec- 
troscopically identified white dwarfs to check kinematically 



distinct spectroscopic subgroups and test different scenar- 
ios of white dwarf production channels. However, a major 
disadvantage of this subset of the white dwarf population is 
that the three components of the velocity are derived only 
from the tangential velocity, since the determination of ra- 
dial velocities for white dwarfs is not an easy task, especially 
for very cool ones. Obviously it would be better to have the 
complete description of the space motions of this sample, 
but it is nonetheless true that we already have two-thirds 
of the motion available for comparison with the simulated 
samples and that the latter samples can account for this 
observational bias. 

The sample of Sion et al. (1988) consists of 626 stars 
with known distances and tangential velocities (of which 
421 white dwarfs belong to the spectral type DA and 205 
stars belong to other spectral types). In this proper mo- 
tion sample there are 523 white dwarfs for which masses, 
radii and effective temperatures could be derived — see Sion 
et al. (1988) for the computational details — of which 372 
have masses larger than 0.5 Mq and, therefore, are expected 
to have carbon-oxygen cores. Of this latter group of white 
dwarfs there are 305 with spectral type DA and 67 belong 
to other spectral types. For this particular sample of white 
dwarfs cooling ages were derived using the cooling sequences 
of Salaris et al. (1997) and, given a relationship between the 
initial mass on the main sequence and the final mass of the 
white dwarf (Iben & Laughlin 1989), main sequence life- 
times (Iben & Laughlin 1989) were also assigned, and the 
birth time of their progenitors was computed. However, the 
errors in the determination of the mass of the progenitor can 
produce large errors in the determination of the total age of 
low mass white dwarfs. For instance, for a typical 0.6 Mq 
white dwarf an error in the determination of its mass of 
0.05 Mq leads to an error in its cooling age of ~ 0.3 Gyr at 
log(L/L Q ) = -2.0 and of ~ 0.8 Gyr at log(L/i Q ) = -4.0, 
whereas the error in the determination of its main sequence 
lifetime is of ~ 2 Gyr. Thus, the mass dependence of the 
cooling sequences is relatively small, whereas the mass de- 
pendence of the main sequence lifetimes is very strong. Fi- 
nally, it could be argued that since this sample includes both 
DA and non-DA white dwarfs, appropiate cooling sequences 
should be used for each spectral type. However, the errors 
introduced by using unappropiate cooling sequences (that is 
cooling sequences for He-dominated white dwarf envelopes) 
in the calculation of the cooling times of DA white dwarfs 
are small when compared to the errors introduced in dating 
white dwarfs by poor mass estimates. Therefore, the tempo- 
ral characteristics of the white dwarf population from them 
derived should be viewed with some caution. Note as well 
that there is not any guaranty that the sample of Sion et 
al. (1988) is representative of the whole population of white 
dwarfs, since it is by no means complete, and therefore some 
cautions are required when drawing conclusions. To be more 
precise, the sample of Sion et al. (1988) has very few low lu- 
minosity white dwarfs. In fact, this sample contains only 
twelve white dwarfs belonging to the low luminosity sample 
of Liebert et al. (1988), of which only four have mass, tan- 
gential velocity, and effective temperature determinations. 
Therefore, we have added to this sample — hereinafter "ob- 
servational" sample — three additional white dwarfs of the 
sample of Liebert et al. (1988) for which a mass estimate 
could be found (Dfaz-Pinto et al. 1994). Nevertheless, this 
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Figure 2. Tangential (upper left panel) and radial (lower left panel) velocity distributions for the original sample and the corresponding 
distributions for the restricted sample (upper and lower right panels, respectively). Also shown as dotted lines are the tangential and 
radial velocity distributions of the restricted sample with a looser restriction in proper motions (see text for details). 



sample provides a unique opportunity to test the results ob- 
tained from a simulated white dwarf sample. 

3.1 The overall kinematical properties of the 
samples 

First we compare the overall kinematical properties of the 
white dwarf simulated samples with those of the observa- 
tional sample, regardless of the birth time of their progen- 
itors. In Figure 2 we show the distributions of the tangen- 
tial and radial velocities for both the original and the re- 
stricted sample. The tangential velocity distribution of the 
original sample is shown in upper left panel of figure 2 and 
the tangential velocity distribution of the restricted sample 
is shown as a solid line (left scale) in the upper right panel. 
The restricted sample, which is kinematically selected, has 
a smaller tangential velocity dispersion (<7tan ~ 80 km s -1 ) 
than the original sample (<7tan ~ 100 km s _1 ). Here we have 
defined for operational purposes only the dispersions to be 
as the full width at half maximum of the distributions. More- 
over, both samples are peaked at different tangential veloc- 
ities: at Kan ~ 45 km s _1 for the original sample and at 
Kan ~ 65 km s _1 for the restricted sample, showing clearly 
that the restricted sample is biased towards larger tangen- 
tial velocities, as it should be for a proper motion selected 
sample. In fact the most probable tangential velocity of the 



restricted sample is almost one third larger than that of the 
total sample of white dwarfs. This kinematical bias is clearly 
seen as well in the behavior of the distribution at low tan- 
gential velocities where the restricted sample shows a deficit 
of low velocity stars, as expected from a kinematically se- 
lected sample. Note as well the existence of an extended tail 
at high tangential velocities, indicating the presence of high 
proper motion white dwarfs. Of course, all these effects are 
simply due to the selection criteria and, in particular to the 
assumed restriction in proper motion. The distribution of 
radial velocities of the original sample is shown in the lower 
left panel of figure 2 and the radial velocity distribution of 
the restricted sample is shown as a solid line (left scale) of 
the lower right panel. Both distributions have similar dis- 
persions (<j ra d ~ 90 km s _1 ) and both are well behaved and 
centered at Kad = 0, as it should be since there is not any 
constrain on the radial velocities of the restricted sample. 

In Figure 3 the tangential velocity distribution of the 
observational sample is shown. By comparing the tangential 
velocity distributions of figure 2 (upper panels) and figure 3 
we can assure that the observational sample does not have 
a clear kinematical bias since it does not show a clear deficit 
of low tangential velocity white dwarfs — the ratio between 
the height of the peak and the height of the lowest velocity 
bin is the same for both the original sample and the ob- 
servational sample: roughly 1/4 — and does not have an 
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Figure 3. Tangential velocity distribution of the sample of Sion 
ct al. (1988). 

extended tail at high tangential velocities as the restricted 
sample does. Moreover, the observational sample peaks at 
Kan ~ 40 km s _1 , whereas the original sample (which is 
not kinematically selected) peaks at a very similar tangen- 
tial velocity (Kan ~ 45 km s _1 ). However the tangential 
velocity dispersion (at an ~ 60 km s 1 ) of the observational 
sample is roughly one-third smaller than that of the original 
sample. This might be due to the absence of low luminos- 
ity white dwarfs in the observational sample. Notice that 
intrisically dim white dwarfs are selected on the basis of a 
large proper motion and, therefore, are expected to have, 
on average, larger tangential velocities, thus increasing the 
velocity dispersion. To check this assumption we have run 
our Monte Carlo code with a looser restriction on proper 
motions (/u > 0.08" yr _1 ). The result is shown in the upper 
right panel of figure 2 as a dotted line (right scale) . Although 
the number of selected white dwarfs increases from ~ 85 to 
almost 250 the tangential velocity dispersion decreases from 
o" t an ~ 80 km s _1 to <7tan ~ 60 km s _1 in good agree- 
ment with the tangential velocity dispersion of the obser- 
vational sample. A final test can be performed by imposing 
a tighter restriction on visual magnitudes (my < 15.5 mag ). 
The resulting sample is now smaller — 58 white dwarfs — as 
should be expected, whereas the tangential velocity disper- 
sion decreases to a ta ,n ~ 40 km s _1 and the most probable 
tangential velocity remains almost unchanged (Kan ~ 40 
km s _1 ). On the other hand the radial velocity distribu- 
tion — dashed line and right scale in the lower right panel 
of figure 2 — is nearly indistinguishable from the previous 
sample, selected with a tighter restriction. Nevertheless, the 
differences between the observational sample and the sim- 
ulated samples could be considered as minor. Therefore we 
conclude that the simulated population of white dwarfs is 
fairly representative of the real population of white dwarfs. 

3.2 The temporal behavior of the samples 

Up to this moment we have compared the global kinemat- 
ical characteristics of the simulated samples with those of 
the sample of Sion et al. (1988), but one of the major ad- 
vantages of this latter sample is that all the mass determi- 



100 I 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 r 




10 8 6 4 2 



ttoui (Gyr) 

Figure 4. Time distributions of the restricted sample (non- 
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of the corresponding bin. 



nations have been obtained using the same procedure, and 
consequently, in this sense, the sample is relatively homo- 
geneous. Therefore, we can tentatively obtain the temporal 
variations of the kinematical properties as a function of the 
birth time of the progenitors of the white dwarfs belonging 
to the observational sample, and compare them with those 
of the simulated samples. 

To this regard, in Figure 4 we show the histograms of 
the distribution of the birth times of white dwarfs belonging 
to the observational sample (shaded histogram) and to the 
restricted sample (non-shaded histogram). The number of 
objects in each time bin is also shown on top of each bin of 
the histogram. Time runs backwards and, therefore, old ob- 
jects are located at the left of the diagrams, whereas young 
objects contribute to the time bins of the right part of the 
diagrams. It is important to realize that old bins may in- 
clude kinematical data coming from either bright, low-mass 
white dwarfs, or dim, massive white dwarfs. The time bins 
have been chosen in such a way that the distribution of white 
dwarfs in the observational sample is efficiently binned. The 
first bin in time corresponds to objects older than 7 Gyr 
and has only 5 objects, most of them corresponding to in- 
trinsically faint objects belonging to the sample of Liebert 
ct al. (1988). The last bin corresponds to objects younger 
than 1 Gyr. The remaining three bins are equally spaced in 
time and correspond to white dwarf progenitors with ages 
running from 1 to 7 Gyr in 2 Gyr intervals. All the bins have 
been centered at the average age of the objects belonging to 
them (~ 7.7, 5.9, 3.7, 1.9 and 0.5 Gyr, respectively). 

Since the youngest time bin corresponds to intrinsically 
bright white dwarfs it is expected that this time bin is rea- 
sonably complete in the observational sample. Therefore we 
have chosen the total number of stars in the simulated sam- 
ples in such a way that the restricted sample has a number 
of objects in the youngest time bin comparable with that 
of the observational sample. Since there is no clear restric- 
tion in the ages of white dwarfs belonging to the restricted 
sample, the statistical reliability of the remaining time bins 
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Figure 5. Components of the tangential velocity as a function of the birth time of the white dwarf progenitors for the restricted (left-hand 
panels) and the observational sample (right-hand panels), see text for details. 



of the observational sample can be readily assessed. Note 
the huge difference in the number of white dwarfs between 
the observational sample and the restricted sample for the 
oldest time bins of figure 4. Clearly, the completeness of 
the observational sample decreases dramatically as the birth 
time increases. The percentage of missing white dwarfs (rj) 
in the observational sample as a function of the birth time 
of their corresponding progenitors is also shown in figure 
4 as a solid line, assuming that the youngest time bin of 
the restricted sample is complete. We have considered the 
observational sample to provide reasonable estimates of the 
temporal variations of the velocity when one third of the 
expected number of white dwarfs is present in the corre- 
sponding time bin. This roughly corresponds to birth times 
smaller than 3.7 Gyr. Therefore, the only time bins that we 



are going to consider statistically significant are the youngest 
three bins. 

In Figure 5 we show as solid lines the temporal variation 
of the components of the tangential velocity as a function of 
the total age (white dwarf cooling age plus main sequence 
lifetime of the corresponding parent star) of white dwarfs 
belonging to the restricted sample (left-hand panels), and 
the same quantities for white dwarfs belonging to the ob- 
servational sample (right-hand panels). Although we have 
not considered the data for ttotai > 3.7 Gyr to be reliable 
due the incompleteness of the observational sample, we also 
show the temporal variations of all the three components of 
the tangential velocity for these times as dotted lines for the 
sake of completeness. The thinner vertical line corresponds 
to ttotai = 3.7 Gyr. As can be seen in this figure, the general 
trend for young objects is very similar for both samples. In 
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Table 1. Average values of the three components of the tangential 
velocity and their corresponding dispersions (both in km s _1 ) for 
several choices of Zf (in kpc). 
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particular both the restricted and the observational sample 
have negative velocities across the galactic plane with veloc- 
ities of W ~ —10km s _1 , both samples lag behind the sun 
with similar velocities of V ~ —25 km s _1 and —20 km s _1 , 
respectively, and both samples have positive radial veloci- 
ties of roughly U ~ 20 km s _1 . Finally, old objects in both 
samples lag behind the sun (middle panels) being the lag ve- 
locity comparable for both samples: V ~ —20 km s _1 and 
—25 km s _1 , respectively. Moreover, we have computed the 
time-averaged values of the velocities shown in figure 5 and 
we have found (U) ~ 10 and 12 km s~\ (V) ~ -28 and 
—23 km s _1 , and (W) ~ —8 and —7 km s _1 , respectively. 
We have computed as well the time-averaged values of the 
velocity dispersions for the restricted and the observational 
samples: (<ttj) ~ 41 and 42 km s _1 , (cry) ~ 27 and 30 
km s _1 , and (aw) ~ 25 and 25 km s _1 , and we have found 
that they are also in good agreement. 

As already noted in §2, the most important ingredient 
needed to fit adequately the kinematics of white dwarfs is the 
exact shape of the scale height law. In fact, the luminosity 
function (see section 4 below) is only sensitive to the ratio of 
the initial to final scale heights (zi/zf) of the disk and to the 
time-scale of disk formation (r h ) but not to the exact value of 
say Zf . However when the kinematics of the sample are con- 
sidered the reverse is true. That is the kinematics of the sim- 
ulated samples are very sensitive to the exact value adopted 
for the final scale height. This is clearly illustrated in Table 
1, where the time-averaged values for the three components 
of the tangential velocity and the tangential velocity disper- 
sions are shown for several choices of the final scale height, 
but keeping constant the above mentioned ratio. As it can 
be seen there, the time-averaged radial component of the 
tangential velocity, (U), and the time-averaged perpendicu- 
lar component of the tangential velocity, (W), are not very 
sensitive to the choice of Zf, whereas the time-averaged lag 
velocity is very sensitive to its choice. Regarding the velocity 
dispersions all three components are sensitive. We have cho- 
sen the value of zt which best fits the average values of the 
observed sample. In order to produce the results of figures 4 
and 5 a value of 500 pc was adopted for Zf, which is typical 
of a thick disk population. It is important to point out here 
that increasing (decreasing) Z{ by a factor of two without 
keeping constant the ratio (zi/zf) doubles (halves) aw for 
objects in the youngest time bin, which is the most reliable 
one, thus making incompatible the simulated and the obser- 
vational samples. Similarly, increasing r h by a factor of two 
changes dramatically the behaviour of the lag velocity since 
it changes the value of V for objects in the youngest time 
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Figure 6. Overall kincmatical properties of the observational 
sample. 



bin from ~ —20 to ~ —10 km s _1 . We conclude that the 
proposed scale height law is not in conflict with the observed 
kinematics of the white dwarf population. 



3.3 A final remark on the reliability of the 
samples 

Finally, it is interesting to compare the results of a kinemati- 
cal analysis of the observational sample with the predictions 
obtained from main sequence star counts as given by equa- 
tion (2). In Figure 6 we show the correlations between the 
V component of the tangential velocity (top panel), the ra- 
tio ov/crv (middle panel), and the ratio <7w/cu (bottom 
panel) as a function of the radial velocity dispersion cttj ob- 
tained from the data of main sequence stars compiled by 
Sion et al. (1988). Except for the last bin the agreement be- 
tween the data obtained from the white dwarf sample and 
the data obtained from main sequence stars is fairly good 
(ffy/^u ~ °-3. °w/^u ~ 0.5 and V/ol ~ -20). However, 
it should be taken into account that the data coming from 
the last bin is obtained, as previously mentioned, with only 
five stars. Moreover, all these objects belong to the low lu- 
minosity sample of Liebert et al. (1988), which is strongly 
biased towards large tangential velocities and, besides, sys- 
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Table 2. Average values of the three components of the tangential 
velocity and their corresponding dispersions (both in km s _1 ) for 
the Monte Carlo simulation, the observational sample, and the 
Edvardsson et al. (1993) sample. 



tematic errors affecting either mass and radius determina- 
tions or luminosity determinations (the bolometric correc- 
tions adopted in the latter work are highly uncertain) can 
mask the true behavior of the sample. 

A final test of the validity of the assumptions adopted 
in this paper to derive the simulated populations can be 
performed by comparing the results of this section with the 
kinematical analysis of a sample of main sequence F and G 
stars (Edvardsson et al. 1993). These authors measured dis- 
tances, proper motions and radial velocities (among other 
data) for a sample of 189 F and G stars. They also assigned 
individual ages for all the stars in the sample from fits in the 
T c ff —log g plane. The same sample has been re-analyzed very 
recently by Ng & Bertelli (1998), using distances based on 
Hipparcos parallaxes and improved isochrones. We refer the 
reader to the latter work for a detailed analysis of the er- 
rors and uncertainties involved in dating individual objects. 
Although an analysis similar to that performed in §3.2 can 
be done, for the sake of conciseness we will only refer here 
to the average values of the three components of the tan- 
gential velocity and its corresponding dispersions. For this 
purpose in Table 2 we show the averaged values of the three 
components of the tangential velocity and their correspond- 
ing dispersions for the restricted sample of our Monte Carlo 
simulation, labelled MC, the observational sample, labeled 
WD, and the three components of the velocity and their 
dispersions for the Edvardsson et al. (1993) sample, labeled 
E93. As already discussed in §3.2 the agreement between 
the Monte Carlo simulation and the observational sample 
is fairly good. The comparison of both samples with the 
sample of Edvardsson et al. (1993) reveals that the agree- 
ment between the average values of the three samples is 
remarkably good, even if the dating procedure for individ- 
ual objects is very different in both observational samples. 
The same holds for the averaged values of the three velocity 
dispersions. We conclude that our equation (2) represents 
fairly well the kinematical properties of the observed white 
dwarf population. 



4 THE WHITE DWARF LUMINOSITY 
FUNCTION 

4.1 The spatial distribution and completeness of 
the simulated white dwarf population 

The 1/Vmax method (Schmidt 1968, Felten 1976), when ap- 
plied to our simulated white dwarf population, should pro- 
vide us with an unbiased estimator of its luminosity func- 
tion, pressumed completeness of the simulated samples in 
both proper motion and apparent magnitude, and provided 




Z (kpc) 

Figure 7. Histogram of the z distribution for both the original 
sample (right scale) and the restricted sample (left scale). 



that the spatial distribution of white dwarfs is homogeneous. 
Strictly speaking this means that the maximum distance at 
which we find an object belonging to the sample is indepen- 
dent of the direction. In our case this is clearly not true — 
and, most probably, for a real sample this would certainly 
be the case as well — since we have derived the simulated 
samples assuming an exponential density profile across the 
galactic plane. Since the scale height law exponentially de- 
creases with time (see §2) it is difficult to say "a priori" 
which is the final spatial configuration of the simulated white 
dwarf samples introduced in the previous sections. 

In the histogram of Figure 7 we show the logarithmic 
distribution of the number of white dwarfs as a function of 
the absolute value of the z coordinate for both the original 
sample (right scale) and the restricted sample (left scale). 
Clearly, both distributions correspond to exponential disk 
profiles with different scale heights. Also shown in figure 
7 are the best fits to these distributions. The corresponding 
scale heights from them derived are ~ 1.3 kpc for the original 
sample, which is typical of a thick disk population, and con- 
siderably smaller « 129 pc for the restricted sample which 
can be considered typical of a thin disk population. This is 
not an evident result since, as has been explained in §2, the 
simulated populations take naturally into account the fact 
that old objects are distributed over larger volumes (that is, 
with larger scale heights and therefore with larger velocity 
dispersion perpendicular to the plane of the galaxy) than 
young ones. Therefore, one could expect that the final spa- 
tial distribution of the restricted white dwarf population — 
which is kinematically selected — should reflect properties 
of an intermediate thin-thick disk population and certainly 
this is not the case. Obviously, since there is not any re- 
striction in the distances (within the local column) at which 
a white dwarf belonging to the original sample can be ob- 
served the expected final scale height for this sample should 
be much larger, in good agreement with the simulations. 
Regarding the restricted sample, our results clearly indicate 
that we are selecting for this sample white dwarfs lying very 
close to the galactic plane. Moreover, if we change by a fac- 
tor of two Zf as explained in §3, the final scale height of 
the restricted sample does not change appreciably and, on 
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Figure 8. Cumulative histograms of apparent magnitude and proper motion for both the original and the restricted sample. See text 
for details. 



the other hand, the dispersion of velocities perpendicular to 
the galactic plane does not agree with its observed value. 
Therefore, the final scale height of the restricted sample is 
clearly dominated by the selection criteria. It is important 
to realize that this scale height, taken at face value, is not 
negligible at all when compared with the value of the max- 
imum distance at which a parallax is likely to be measured 
with relatively good accuracy — which is typically 200 pc — 
and which imposes an additional selection criterion (see §2) 
for white dwarfs belonging to the restricted sample, which 
are the white dwarfs which are going to be used in the pro- 
cess of determination of the white dwarf luminosity function. 
Therefore, the l/Vmax method must be generalized to take 
into account a space-density gradient. For this reason we 
have used the density law of figure 7 to define a new density- 
weighted volume element dV' = p(z) dV (Felten 1976; Avni 
& Bahcall 1980; Tinney, Reid & Mould 1993), being p(z) the 
density law derived from figure 7. This new, corrected, esti- 
mator provides a more accurate determination of the white 
dwarf luminosity function and, ultimately, a more realistic 
value of the space density of white dwarfs. All in all, for 
reasonable choices of a scale height law its effects on the 
derived white dwarf luminosity function in principle cannot 
be considered negligible. 

The second, and probably more important issue, is the 
completeness of the samples used to build the white dwarf 
luminosity function. This is a central issue since the 1/Knax 



method assumes completeness of the samples. The reader 
should keep in mind that the original sample is complete 
by construction, since it consists of all white dwarfs gener- 
ated by the Monte Carlo code, regardless of their distance, 
proper motion, apparent magnitude and tangential veloc- 
ity; whereas the restricted sample is built with white dwarfs 
culled from the original sample according to a set of selec- 
tion criteria and, therefore, its completeness remains to be 
assessed. 

In Figure 8 we explore the completeness of the simu- 
lated samples. For this purpose, the cumulative star counts 
of white dwarfs with apparent magnitude smaller than rnv 
for the original sample are shown in the top left panel of fig- 
ure 8, whereas the corresponding diagram for the restricted 
sample is shown in the top right panel. Also shown in fig- 
ure 8 are the cumulative star counts of white dwarfs with 
proper motions larger than fi belonging to the original sam- 
ple (bottom left panel) and to the restricted sample (bottom 
right panel) . For a complete sample distributed according to 
a homogenous spatial density, the logarithm of the cumula- 
tive star counts of white dwarfs with apparent magnitude 
smaller than mv are proportional to mv with a slope of 0.6 
(see, for instance, Mihalas & Binney 1981). We also show in 
the top panels of figure 8 a straight line with such a slope. It 
is evident from the previous discussion that our samples are 
not, by any means, distributed homogenously. Note as well 
that the effects of a scale height law are tangled in the stan- 
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dard test of completeness of the samples. Nevertheless, the 
effects of a scale height law can be disentangled since they 
should be quite apparent in the cumulative star counts dia- 
gram of the original sample, which is complete. A look at the 
top left panel of figure 8 reveals that the effects of the scale 
height law are evident for surveys with limiting magnitude 
mv > 19 mag . Therefore, we can now assess the completeness 
in apparent magnitude of the restricted sample, since the 
turn-off for this sample (see top right panel of figure 8) oc- 
curs at my ~ f7 mag . Consequently, the effects of the scale 
height law can be completely ruled out, and this value can 
be considered as a safe limit for which the restricted sample 
is complete in apparent magnitude. 

The completeness of the restricted sample in proper mo- 
tion can be assessed in a similar way. Again, the assumption 
of an homogenous and complete sample in proper motion 
leads to the conclusion that the logarithm of the cumulative 
star counts of white dwarfs with proper motion larger than 
H should be proportional to /i with a slope of —3 (see, for in- 
stance, Oswalt & Smith 1995, and Wood & Oswalt 1998). A 
look at the bottom left panel of figure 8 reveals that for the 
original sample this is not by far the case. In other words, 
since this particular sample is complete by construction the 
hypothesis of an homogenous distribution of proper motions 
must be dropped. This is again one, and probably the most 
important, of the effects associated with a scale height law 
since the kinematics of the samples are highly sensitive to 
the choice of the scale height law (see §2, equation (2) and 
table 1). It is important to realize that the effects of a scale 
height law are more prominent in proper motion than in the 
spatial distribution and this can be directly checked for a 
real sample, thus providing a direct probe of the history of 
the star formation rate per unit volume. Finally, in the lower 
right panel of figure 8 the cumulative star counts in proper 
motion of white dwarfs belonging to the restricted sample 
are shown. As expected, the effects of a scale height law are 
in this case negligible since we are culling white dwarfs with 
high proper motion for which the original sample is reason- 
ably complete (see the lower left panel of figure 8). The exact 
value of the turn-off is in this case /i ~ 0.3" yr _1 , in close 
agreement with the results of Wood & Oswalt (1998). 

It is quite clear from the previous discussions that one 
of the ingredients that has proven to be essential in the de- 
termination of the white dwarf luminosity function is the 
adopted scale height law. In principle one should expect two 
kinds of competing trends. On the one hand, the effects of 
the scale height law should be more dramatic for old ob- 
jects, because old objects have a larger velocity dispersion 
(although the effects of a spatial inhomogeneity should be, 
as well, less apparent) and the tail of the white dwarf lumi- 
nosity function is populated predominantly by this kind of 
white dwarfs (intrinsically dim, high proper motion objects). 
Therefore, one should expect that the cut-off in the white 
dwarf luminosity function is influenced either by the spatial 
distribution of white dwarfs or by their velocity distribution 
or by a combination of both. On the other hand, objects pop- 
ulating the tail of the luminosity function are intrinsically 
dim objects and, therefore, in order to be selected for the re- 
stricted sample they must be close neighbors. This, in turn, 
implies that the average distance at which we are looking 
for white dwarfs is small and, consequently, the effects of a 
scale height law should be less apparent. The reverse is true 
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Figure 9. Average properties of the restricted sample as a func- 
tion of the luminosity. Bottom panel: average z coordinate of 
white dwarfs of the restricted sample; middle panel: average tan- 
gential velocity for these white dwarfs, the observational data 
has been taken from Liebert et al. (1988); and top panel: average 
proper motion of these objects. 



at moderately high luminosities. Therefore it is interesting 
to see which are the dominant effects as a function of the 
luminosity. For this purpose, in Figure 9 we show several av- 
erage properties of white dwarfs belonging to the restricted 
sample as a function of their luminosity for a typical Monte 
Carlo simulation. 

In the top panel of figure 9 we show the average dis- 
tance to the galactic plane of white dwarfs belonging to the 
restricted sample. As it can be seen, the average distance to 
the galactic plane of intrinsically bright white dwarfs can be 
as high as 100 pc, which is a sizeable fraction of the derived 
scale heights of the white dwarf samples. Consequently we 
expect that the effects of an inhomogeneous spatial distri- 
bution should be very prominent at high luminosities. Con- 
versely, the average distance to the galactic plane for white 
dwarfs near the observed cut-off in the white dwarf lumi- 
nosity function is only ~ 10 pc. Therefore as the luminosity 
decreases we are probing smaller volumes and the effects 
of an inhomogenous spatial distribution at low luminosities 
are expected to be, from this point of view, small (but see, 
however, the discussion in §4.4) . 

In the middle panel of figure 9 the average tangential ve- 
locity of objects belonging to the restricted sample is shown 
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as a function of the luminosity. The observational data is 
shown as solid circles and has been obtained from Liebert 
et al. (1988). Their adopted restriction in proper motion 
fM> = 0.80" yr _1 is significantly larger than the one we adopt 
here /io = 0.16" yr~\ which is consistent with the cut-off in 
proper motion adopted by Oswalt et al. (1996). Therefore we 
expect a smaller average tangential velocity. The agreement 
is fairly good since, given the ratio of proper motion cut- 
offs, the average tangential velocity of our restricted sample 
should be roughly a 20% smaller: a closer look at the middle 
panel of figure 9 shows that the average tangential velocity 
reported by Liebert et al. (1988) is ~ 120 km s _1 , whereas 
we obtain ~ 90 km s _1 . These figures reinforce the gen- 
eral idea that our simulations are fully consistent with the 
observed kinematics of the white dwarf population. 

Finally, in the bottom panel of figure 9 the average 
proper motion distribution of those stars belonging to the 
restricted sample is shown as a function of the luminosity. 
As it can be seen there, low luminosity white dwarfs belong- 
ing to the restricted sample have, on average, large proper 
motions. As expected from the discussion of the two previ- 
ous panels, the distribution of proper motions is smoothly 
increasing for luminosities in excess of ~ 10~ 3 Lq: since 
the average value of the tangential velocity remains approx- 
imately constant, and we are selecting objects with smaller 
average distances, the net result is an increase in the average 
proper motion. Moreover, white dwarfs belonging the low 
luminosity portion of the white dwarf luminosity function 
are preferentially culled from the original sample because 
of their high proper motion. That is the same to say that 
the selection criterion is primarily the proper motion one 
and that the criterion on apparent magnitude has little to 
do for these luminosities, in agreement with the results of 
Wood & Oswalt (1988). As a final consequence the effects 
of an inhomogenous distribution in proper motions will be 
more evident at high luminosities, where the average proper 
motion is smaller (see the discussion of the lower left panel 
of figure 8). All in all, the effects of the inhomogeneities in 
both proper motion and z will be more prominent at high 
luminosities, where the observational luminosity function al- 
ready takes into account these effects (Fleming et al. 1986). 

4.2 The Monte Carlo simulated white dwarf 
luminosity functions 

In Figure 10 we show a set of panels containing the 
white dwarf luminosity functions obtained from ten different 
Monte Carlo simulations. That is the same to say that ten 
different initial seeds were chosen for the random number 
generator and, consequently, ten independent realizations 
of the white dwarf luminosity function were computed (in 
fact, we have computed twenty independent realizations, of 
which only ten are shown in figure 10). The adopted age of 
the disk was fdisk = 13 Gyr and the set of restrictions used 
to build the sample is that of §2, which is the same set used 
by Oswalt et al. (1996) to derive their observational white 
dwarf luminosity function. The simulated white dwarf lumi- 
nosity functions were computed using a generalized 1/V m ax 
method (Felten 1976; Tinney et al. 1993; Qin & Xie 1997) 
which takes into account the effects of the scale height. The 
error bars of each bin were computed according to Liebert 
ct al. (1988): the contribution of each star to the total error 
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Figure 10. Panel showing different realizations of the simulated 
white dwarf luminosity function — filled squares and solid lines 
— compared to the observational luminosity function of Oswalt 
et al. (1996) — filled circles and dotted line. 



budget in its luminosity bin is conservatively estimated to be 
the same amount that contributes to the resulting density; 
the partial contributions of each star in the bin are squared 
and then added, the final error is the square root of this 
value. The resulting white dwarf luminosity functions are 
plotted as solid squares; a solid line linking each one of their 
points is also shown as a visual help. Also plotted in each 
one of the panels is the observational white dwarf luminosity 
function of Oswalt et al. (1996), which is shown as solid cir- 
cles linked by a dotted line. For each realization of the white 
dwarf luminosity function the value obtained for (V/Vmax) 
is also shown in the upper left corner of the corresponding 
panel. Finally, and for sake of completeness, the total num- 
ber of objects in the restricted sample, JVwd, of the different 
realizations of the white dwarf population, and the distribu- 
tion of objects, Ni, in each luminosity bin, i, of the Monte 
Carlo simulated white dwarf luminosity functions are shown 
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Table 3. Total number of white dwarfs, A?wd> an d white dwarfs in each bin, AT;, for each of the twenty realizations of the simulated 
white dwarf luminosity functions. 



Table 4. 
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1 
2 
3 
4 
5 
6 
7 
8 
9 
10 
11 
12 
13 
14 
15 
16 
17 
18 
19 



0.24 



0.80 
0.46 



0.21 
0.09 
0.40 



0.41 
0.93 
0.69 
0.63 



0.37 
0.26 
0.81 
0.86 
0.65 



0.05 
0.31 
0.78 
0.13 
0.35 
0.77 



0.03 
0.59 
0.22 
0.01 
0.16 
0.01 
0.11 



0.16 
0.48 
0.73 
0.81 
0.80 
0.92 
0.68 
0.12 



0.24 
0.79 
0.69 
0.20 
0.90 
0.32 
0.52 
0.67 
0.58 



0.33 
0.48 
0.88 
0.88 
0.92 
0.89 
0.47 
0.11 
0.99 
0.54 



0.32 
0.30 
0.90 
0.19 
0.38 
0.86 
0.89 
0.23 
0.77 
0.69 
0.56 



0.04 
0.37 
0.44 
0.24 
0.26 
0.59 
0.80 
0.20 
0.82 
0.58 
0.41 
0.79 



0.70 
0.36 
0.89 
0.77 
0.68 
0.94 
0.55 
0.03 
0.77 
0.54 
0.72 
0.84 
0.62 



0.86 
0.47 
0.85 
0.60 
0.65 
0.85 
0.34 
0.06 
0.68 
0.57 
0.67 
0.67 
0.51 
0.79 



0.17 
0.29 
0.67 
0.14 
0.39 
0.44 
0.27 
0.46 
0.77 
0.45 
0.85 
0.77 
0.38 
0.24 
0.47 



0.06 
0.11 
0.56 
0.57 
0.43 
0.89 
0.61 
0.04 
0.83 
0.51 
0.81 
0.77 
0.60 
0.38 
0.65 
0.78 



0.75 
0.67 
0.81 
0.36 
0.90 
0.32 
0.11 
0.31 
0.53 
0.74 
0.87 
0.29 
0.15 
0.40 
0.82 
0.79 
0.38 



0.02 
0.19 
0.16 
0.37 
0.54 
0.83 
0.42 
0.01 
0.59 
0.07 
0.46 
0.15 
0.10 
0.25 
0.16 
0.04 
0.28 
0.04 



0.34 
0.86 
0.89 
0.26 
0.84 
0.49 
0.50 
0.76 
0.92 
0.81 
0.95 
0.71 
0.63 
0.56 
0.61 
0.96 
0.50 
0.93 
0.12 



in Table 3. The total number of white dwarfs belonging to 
the restricted sample is roughly 200, which is the typical 
size of the samples used to build the currently available ob- 
servational luminosity functions. This number is important 
since the assigned error bars are strongly dependent on the 
number of objects in each luminosity bin. 

It is important to notice the overall excellent agreement 
between the simulated data and the observational luminos- 
ity function. However, there are several points that deserve 
further comments. The first one is that the simulated white 
dwarf luminosity functions are systematically larger than 
the observational luminosity functions for luminosities in ex- 
cess of \og(L/ Lq) = —2.0. This behavior reflects the effects 
of the spatial inhomogeneity ot the simulated white dwarf 



samples. It is important to realize that the hot portion of 
the white dwarf luminosity function of Oswalt et al. (1996) 
has been derived without taking into account the effects 
of a scale height, in contrast with the procedure adopted 
by Fleming et al. (1986), where those effects were properly 
taken into account. When one compares the luminosity func- 
tions obtained in this section with that of Fleming et al. 
(1986) the agreement is excellent. Also of interest is the fact 
that the hot portion of the white dwarf luminosity func- 
tion varies quite considerably for the different realizations. 
The reason for this behavior is that at high luminosities 
the evolution is dominated by neutrino losses and it is fast. 
Therefore, the probability of finding such white dwarfs is 
relatively small and the statistical signifance of those bins is 
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Figure 11. Probability distribution functions for each luminosity bin. 



low. Consequently, the exact shape of the luminosity func- 
tion at \og(L/ Lq) > —3.0 is strongly dependent of the ini- 
tial seed of the pseudo-random number generator. This is 
further confirmed by comparing the second and the third 
column of table 3 where the total number of objects in the 
restricted sample and the number of objects in the first bin 
of the white dwarf luminosity function of each realization of 
the Monte Carlo simulations are shown. As a consequence 
the real error bars that should be assigned to each bin are 
presumably larger than those of figure 10. Moreover, any 
attempt to derive the volumetric star formation rate using 
data from the bins at high luminosities (Noh & Scalo, 1990) 
is based on very weak grounds. It is also important to notice 
that the completeness of the simulated samples as derived 
from the value of (V/Vmax) is relatively large. In fact, for 
a complete and homogeneous sample this value should be 
equal to 0.5; since the simulated sample samples are not 



homogenous the values obtained here can be considered as 
reasonable. 

Finally, it is convenient to point out here that we have 
done a x 2 test of the compatibility of the Monte Carlo sim- 
ulated samples. The results are shown in Table 4, where the 
probability of an independent observer to find the realiza- 
tions compatible is shown for each pair of realizations. As 
it can be seen, this probability can be as low as 0.01, which 
is the same to say that the corresponding luminosity func- 
tions are completely incompatible, even if they have derived 
from the same set of input parameters and selection criteria. 
Obviously, the conclusion is that for a reasonable number of 
objects in the restricted sample, the white dwarf luminosity 
function is dominated by the selection criteria. 
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Table 5. Error bars of the twenty independent realizations of the Monte Carlo simulated white dwarf luminosity functions in each 
luminosity bin and the same quantities for the Bayesian luminosity function (last row). 
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A log(n,2) 
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A log(ns) 


A log(ri6) 


A log(ri7) 


A log(ns) 


A log(ng) 
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1.000 
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0.000 
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0.387 
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0.500 


0.304 


0.205 


0.439 


10 


1.000 


1.000 


0.541 


0.440 


0.496 


0.275 


0.213 


0.195 


0.628 


11 


0.000 


0.763 


0.577 


0.312 


0.258 
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0.381 


12 


1.000 


0.746 


0.441 


0.446 


0.404 
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13 


0.000 
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0.450 


0.362 
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14 


0.000 


0.741 
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0.486 


0.281 


0.344 


0.195 


0.190 


0.425 


15 
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0.887 


0.560 
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0.292 


0.357 
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0.178 
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16 


1.000 


0.890 


0.450 
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0.335 


0.344 
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17 


0.785 
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0.446 
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1.000 


0.537 


0.496 


0.536 


0.349 


0.254 


0.197 


0.161 


0.484 


19 


1.000 
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— oo 
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4.3 A bayesian analysis of the simulated samples 

As previously stated, changing the initial seed of the random 
function generator the Monte Carlo code provides different 
independent realizations of the white dwarf luminosity func- 
tion. All these realizations are "a priori" equally good. Be- 
sides, since the number of objects that is used to compute 
the white dwarf luminosity function is relatively small, large 
deviations are expected, especially at relatively high lumi- 
nosities for which the cooling timescales are short. This, in 
turn, results in very probable underestimates of the asso- 
ciated uncertainties, especially at luminosities larger than 
log(L/I/Q) <; —3.0. Consequently we have used bayesian 
statistical methods (Press 1996) to obtain a realistic esti- 
mation of the errors involved and the most probable value 
of the density of white dwarfs for each luminosity bin. 

The problem can be stated as follows: for a given lumi- 
nosity, L, we want to know the most probable value of the 
white dwarf luminosity function, N, given a set of Ni sim- 
ulations assuming that all simulations are equally good. To 
compute N one must maximize the probability distribution 



(3) 



where Pg and Pb are the probability of being a good and a 
bad simulation, respectively. We can calculate them follow- 
ing closely Press (1996): 



Pg; = exp 



Pbj = exp 



(Nj - N) 2 
2a? 



(Ni - N) 2 
2S 2 



(4) 



where <7 ; is the error bar of each bin of the luminosity func- 
tion and S is a large but finite number characterizing the 



maximum expected deviation in Ni. We recall that the con- 
tribution to the error of each white dwarf is equal to the 
inverse of its maximum volume squared. 

The results are shown in Figure 11, where the prob- 
ability distributions corresponding to each luminosity bin, 
computed with the previous method, are displayed. The log- 
arithm of the luminosity of each bin in solar units is shown in 
the upper right corner of each panel. All the probability dis- 
tributions, except that of the brighter luminosity bin, have 
a Gaussian profile. This is a direct consequence of the poor 
statistical significance of the first bin. In order to produce 
these probability distributions 20 independent realizations 
of the simulated samples were used. This is a reasonable 
number: increasing the total number of simulations does not 
introduces substantial improvements in the statistical sig- 
nificance of the first bin, which is the less significant. From 
these probability distributions a better estimate of the sta- 
tistical noise can be obtained. We have estimated the result- 
ing error bars by assuming a conservative 95% confidence 
level (approximately 2a) . In Table 5 we show the computed 
deviations for each of the twenty realizations of the Monte 
Carlo simulated white dwarf luminosity functions and the 
most probable error bars computed at the 95% confidence 
level. The error bars obtained from a bayesian analysis of 
the twenty Monte Carlo simulations compare favourably, 
roughly speaking, with those of each individual Monte Carlo 
simulation. However, for samples where the total number of 
white dwarfs is smaller than 200 (the simulations presented 
here) the errors for each of the luminosity bins are severe 
underestimates of the real errors, especially at low luminosi- 
ties. 

In Figure 12 the most probable white dwarf luminos- 
ity function — hereinafter bayesian white dwarf luminosity 
function — with its corresponding error bars is shown, ob- 
tained by maximizing the probability distributions of figure 
11. Except for moderately high luminosities — i.e. for lumi- 
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Figure 12. Bayesian luminosity function. 



nosities larger than log(L/Z/Q) = —2.0 — where the effects 
of the spatial inhomogeinities are most obvious the agree- 
ment between the observational luminosity function and the 
bayesian luminosity function is excellent. Moreover, for the 
bayesian white dwarf luminosity we have computed a syn- 
thetic value of (V/Vmax) as an average of the corresponding 
values for each of the twenty realizations with the weights 
given by the probability of each realization obtained from 
the probability distributions of figure 11. We have obtained 
a value of (V/Vmax) = 0.464 which remains close to the 
canonical value of (V/Knax) = 0.5, valid for an homogenous 
and complete sample. 



4.4 The age of the disk 

Perhaps one of the most surprising results of the simulations 
presented here is the age of the disk itself. The value of 13 
Gyr adopted in this paper fits nicely the observational data 
of Oswalt et al. (1996) as can be seen in figure 12. This 
a direct consequence of the adopted scale height law, since 
using the same set of cooling sequences and a conventional 
approach to compute the white dwarf luminosity function 
with a constant volumetric star formation rate, Salaris et 
al. (1997) derived an age for the solar neighborhood of 11 
Gyr when the effect of phase separation upon crystallization 
was taken into account and of 10 Gyr when phase separation 
was neglected. Thus the ultimate reason of the increase in 
the adopted age of the solar neighborhood is not due to 
the details of the adopted cooling sequences. Instead, this 
increase can be easily explained in terms of the model of 
galactic evolution. We recall that the white dwarf luminosity 
function measures the number of white dwarfs per cubic 
parsec and unit bolometric magnitude. Therefore in order 
to evaluate it the volumetric star formation rate is required. 
In our case we can define the effective star formation rate per 
cubic parsec as ip e s{t) ~ ip(t)/H p (t). With the laws adopted 
here for ip(t) and H p (t) it is easy to verify that the effective 
star formation rate only becomes significant after ~ 2 Gyr 
(Isern et al. 1995a,b). 
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Figure 13. Panel showing different simulated luminosity func- 
tions (filled squares) and their corresponding fit using a standard 
method. 



4.5 Statistical uncertainties in the derived age of 
the disk 

The easiest and more straightforward way to assess the sta- 
tistical errors associated with the measurement of the age of 
the solar neighborhood is trying to reproduce the standard 
procedure. That is, we have fitted the position of the "ob- 
servational" cut-off of each of the Monte Carlo realizations 
with a standard method (Hernanz et al. 1994) to compute 
the white dwarf luminosity function using exactly the same 
inputs adopted to simulate the Monte Carlo realizations, ex- 
cept, of course, the age of the disk, which is the only free 
parameter. The results are shown in Figure 13 for ten of 
the twenty realizations. As is usual with real observational 
luminosity functions the theoretical white dwarf luminos- 
ity functions were normalized to the bin with minimum er- 
ror bars. The derived ages of the disk for each one of the 
realizations are shown in the upper left corner of the cor- 
responding panel. As it can be seen, there is a clear bias: 
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Figure 14. Bayesian analysis of the derived age of the disk. 

the derived ages of the disk are systematically larger than 
the input age of the Monte Carlo simulator by about half 
a Gyr. This a direct consequence of the binning procedure, 
since we are grouping white dwarfs belonging to the max- 
imum of the white dwarf luminosity function in the lowest 
luminosity bin, and can be avoided by using the cumulative 
white dwarf luminosity function, which minimizes the effects 
of the binning procedure. 

We have used the bayesian inference techniques de- 
scribed in §4.3 to assign a purely statistical error to our age 
estimates. In order to do this we need to know a formal un- 
certainty for each one of the independent realizations. Since 
the value of {V/Vma,*} is a good measure of the overall quality 
of the sample (despite the fact that the samples are inho- 
mogenous) we have adopted <7i = 2(0.5 — (V/Vmax)) ^disk- 
The correspoding probability distribution is shown in Fig- 
ure 14, which leads to a most probable age of the disk of 
tdisk = 13.5 ± 0.8 Gyr at the 95% confidence level (2er) 
which has to be compared with the adopted input age in 
our Monte Carlo simulation which was taken to be 13 Gyr. 
This uncertainty is in good agreement with the results of 
Wood & Oswalt (1998). If one relaxes the confidence level 
to lcr the associated error bar is ± 0.4 Gyr. It is nonetheless 
important to realize that there is a systematic increase in 
the inferred disk ages of ~ 5%. 



5 CONCLUSIONS 

We have built a Monte Carlo code which has allowed us 
to reproduce quite accurately the process of building the 
white dwarf luminosity function from the stage of object 
selection and data binning. Our Monte Carlo simulation in- 
cludes a model of galactic evolution based on well established 
grounds. The simulated samples obtained with our code re- 
produce very well the overall kinematical properties of the 
observed white dwarf population. We have explored as well 
the temporal evolution of such properties and we have found 
that there is a fair agreement between the observed distri- 
butions and the simulated ones. However, we have shown 
that observed sample has very few old white dwarfs and 
this has constrained our kinematical analysis to the most 



recent 3.7 Gyr, which can be considered as reasonably se- 
cure. Nevertheless, we have shown that the precise shape of 
the temporal distributions encode a wealth of information 
and, therefore, more detailed analysis should be undertaken 
with improved observational samples. We have also extended 
our kinematical analysis by comparing the averaged proper- 
ties of both the observational white dwarf sample and the 
Monte Carlo simulated sample with a sample of old main 
sequence F and G stars. The results of this comparison lead 
to the conclusion that our model of galactic evolution is fully 
compatible with the properties of the observed samples. 

Using our synthetic populations we have assessed the 
completeness of the samples used to derive the white dwarf 
luminosity function and studied their spatial distribution, 
and given a set of selection criteria consistent with the ob- 
servational procedures, we have built several independent 
realizations of the white dwarf luminosity function and com- 
pared them to the observational luminosity function of Os- 
walt et al. (1996). Our results regarding the white dwarf 
luminosity function can be summarized as follows: 

(i) Given the selection criteria adopted by Oswalt et al. 
(1996), our Monte Carlo simulation strongly suggest that the 
observational samples are complete up to 17 mag and that the 
primary selection criterion at low luminosities is the proper 
motion one, in agreement with Wood & Oswalt (1998). 

(ii) The Monte Carlo simulated white dwarf luminosity 
functions present an excellent agreement with the observa- 
tional data. 

(iii) The effects of a scale height law are not negligible 
at all in the spatial distribution of the samples, especially 
at moderately large luminosities were they are more promi- 
nent. A scale height of roughly 130 pc is derived from the 
Monte Carlo simulations for the objects used in building 
the white dwarf luminosity function. However, we have es- 
tablished without any doubt that the effects of a scale height 
law should be more apparent in the cumulative distribution 
of proper motions. Although the effects of the scale height 
law on the tail of the white dwarf luminosity function seem 
to be negligible at first glance, a detailed analysis reveals 
that this inflation effect increases the derived ages of the disk 
by a considerable amount, which can be typically 2 Gyr. 

(iv) By using bayesian inference techniques we have been 
able to establish that the current procedure to assign the ob- 
servational error bars to the white dwarf luminosity function 
is reasonably for a sample of 200 stars. 

(v) Finally, the statistical uncertainty in the age of the 
disk derived from a bayesian analysis is roughly 1 Gyr, in 
agreement with Wood & Oswalt (1998), and we have deter- 
mined that there is a systematic trend due to the binning 
procedure which increases the disk ages inferred from the 
observational luminosity function by roughly a 5%. 

Nevertheless, a good deal of work remains to be done. 
Future improvements may include a more detailed analy- 
sis of the kinematical properties of the sample of old white 
dwarfs. For this purpose, it would be very useful to have 
more reliable observational samples. In our case this means 
not only complete samples but also more accurate mass de- 
terminations. It would be also convenient to analyze the 
three dimensional motion of the white dwarf population but, 
for the moment, this seems to be unavoidable for the faintest 
white dwarfs due to the absence of spectral features. Also of 
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interest is to study the contamination of the input samples 
used in the process of building the white dwarf luminosity 
function with white dwarfs belonging to the galactic halo. 
This, in principle, cannot be discarded since halo members 
are selected on the basis of high proper motion, which is 
the dominant selection criterion at the dim end of the disk 
white dwarf luminosity function. A detailed statistical anal- 
ysis of the cumulative counts of white dwarfs still remains 
to be done, instead of using the differential space density. 
Last but not least, the tests proposed in this paper could be 
applied to a real sample, thus providing us with very useful 
hints about the structure and evolution of our Galaxy. 
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